Full counting statistics in the self-dual interacting resonant level model 
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We present a general technique to obtain the zero temperature full counting statistics of charge 
transfer in interacting impurity models out of equilibrium from time-dependent simulations on a 
lattice. We demonstrate the technique with application to the self-dual interacting resonant level 
model, where very good agreement between numerical simulations using the density matrix renor- 
malization group and those obtained analytically from the thermodynamic Bethe ansatz is found. 
We show from the exact form of counting statistics that the quasiparticles involved in transport 
carry charge 2e in the low bias regime, and e/2 in the high bias regime. 

PACS numbers: 73.63.-b, 72.70.+m, 05.40.Ca, 05.60.Gg 



As nano-devices decrease in size and electron corre- 
lations become more important, the study of out-of- 
equilibrium transport properties on the nanoscale be- 
comes of fundamental importance. Of great worth are 
powerful numerical approaches to interacting many-body 
systems, such as the density matrix renormalization 
group (DMRG) [l[ . While recent work using the DMRG 
has shown remarkable development in calculation of I-V 
Q and noise 0] characteristics of impurity models, the 
grail of such an enterprise is to calculate all moments of 
charge transfer when the system is driven out of equilib- 
rium, known as the full counting statistics (FCS) 

The electron FCS in mesoscopic transport concentrates 
on a charge distribution p n , which is the probability that 
n electrons are transmitted from the left to the right lead 
during the measurement time t m . Knowledge of all prob- 
abilities p n gives full knowledge of the statistics of charge 
transfer. It is usually more convenient to work with the 
cumulant generating function (CGF) of the distribution 
which depends on the counting field \ and is defined as 
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The last equality here is well defined in the classical limit, 
when the charge Q is a random variable associated with 
some stochastic process [6[ . In the quantum case, one has 
to supplement this definition with a prescription for how 
the measurements are time-ordered; it was shown in j4| 
that in terms of a simple model of a spin-i galvanometer 
as a measurement device that the correct prescription is 

Z t J x ) = (Tte 1 ^ Jo™ TV* /o m (2) 

where T t means time-ordered, T t means anti-time or- 
dered, and the average (. . .) is taken over some ini- 
tial (in general non-equilibrium) state. We deal with 
F(x) = — lnZ(x), as this generates the irreducible cu- 



the long-time limit, each of these is proportional to the 
measurement time t m , for example C\ = It m where / is 
the current flowing in the system, and Ci = St m where 
S is the zero-frequency shot noise. 

While the FCS for non-interacting electrons is fairly 
comprehensively understood J4 L 7[ , far less is known in 
the presence of interactions [8|, \9 1 , particularly when the 
interactions drive the system into a strongly correlated 
regime. In this paper we describe a general technique to 
obtain the CGF for quantum impurity models at finite 
bias by means of time-dependent numerical simulation. 

As a concrete example, we take the interacting reso- 
nant level model (IRLM), described by the Hamiltonian 

L/2 

H = -t ( c L c «.m + H - c ) + ( £ o - U)dU (3) 

n=L,R i=0 



Here, c„ i creates a Fermion on the i'th site of the left 
or right lead, while <fi is the creation operator on the 
resonant level. The parameter t = 1 is the hopping pa- 
rameter of the leads, eo is the energy of the resonant 
level, t' n is the hybridization between the resonant level 
and the leads (here, we will assume a symmetric coupling 
t' R = t' L ), and U is the interaction between the resonant 
level and the leads. In the present work, we concentrate 
on the point U = 2t where the model shows a certain 
self-duality @, [l(| ■ The separation of the above Hamilto- 
nian into two leads and the quantum impurity however 
is much more general than this model and the method 
described below may be readily applied to other setups. 

The counting field is added via the substitution 
t' L r R \ — > i'e ±JX / 4 and if we define the resulting Hamilto- 
nian after this substitution as 'H x , it can be shown that 
the CGF defined in Eq. [2] may be rewritten as 
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FIG. 1: [Color online] Time evolution of the imaginary part 
of — dt m F(x)/(xVsD) for various values of bias voltage Vsd, 
system size M, and small values of counting field \- As ex- 
plained in the text, this converges in the long time limit to 
the conductance of the model. The analytic result for conduc- 
tance is also plotted. It is seen that after an initial transient, 
excellent agreement is seen between the numerical simulation 
in the presence of the counting field and the analytic result. 



This is now very convenient for a numerical time evo- 
lution approach, as \^(t m )) = e~' Mxtm |\P(0)) is sim- 
ply the time evolution of the starting state ^(0)) - 
see also Ref. 0- Following Ref. [ll|, we prepare an ini- 
tial state |^(0)) at time r = by applying a poten- 
tial ±Vsd/2 to the left and right leads respectively. We 
then switch off the voltage bias and evolve (still with 
X = 0) until some time to, when the largest transients 
have died away, and we are roughly in a steady state. 
We then switch on \ and perform two time evolutions 
|*± x (r)) =exp(-i(H± x -E )(T-t ))\y(t )), where E 
is the eigen-energy of the initial state. After each time 
step, we look at the overlap between these two states 
(\&_ x (t)|\I> +x (t)), thus producing from Eq. H]Z(x) as a 
function of the measuring time t rn = t — t$. To obtain 
the long time t m — > oo steady state F(x) = t m F(x) from 
the finite time numerical data, we fit to 



dF( X ) 
dt m 
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-7), 
(5) 

where d gj^ is obtained by the derivative of cubic splines 
fitted to F(x), taking care with the complex logarithm 
so as to get a continuous function. Here, the first cosine 
term is the Josephson current that appears from a finite 
size gap Q, while the final term is a decaying transient. 

We now look at the results of applying this method to 
the self-dual IRLM, defined by Hamiltonian ([3]). For the 
time-evolution, we use the full time-dependent extension 
12j of the DMRG where we optimize the target space 
for all time steps simultaneously avoiding the run-away 
errors in the wave function in adaptive time evolution 
schemes 13(. The time stepping is performed by a ma- 



FIG. 2: [Color online] Time evolution of the real part of 
dt m 2F(x)/x 2 f° r various values of Vsd, M and small values 
of x- This gives the zero frequency shot-noise of the model, 
which is known to have a noticeable finite-size correction (fsc) 
3]. Excellent agreement is found between the time evolved 
simulation and the fsc analytic result, which is also plotted. 



trix exponential using the Krylov space representation 
12]. We first perform short time dynamics up to t m = 2 
and successively restart the DMRG procedure to longer 
measuring times up to t m = 25 using time steps of 0.25. 
We keep up to 3100 states per DMRG block leading to 
target space dimensions of up to 8 • 10 6 . We would like to 
emphasize that the general framework is not tied to the 
DMRG and could be applied within other numerical ap- 
proaches, or even to bosonic problem such as the photon 
dynamics in wave guiding structures |14j . 

We begin the discussion by looking at the overlap 
for small values of x, where one can write F(x) ~ 
—iCiX + C2X 2 /2- Hence the imaginary part at x < 1 
gives the current. We plot this in Fig. [T] for various val- 
ues of bias voltage Vsd and counting fields %, where we 
see relatively weak transient effects, and a clear conver- 
gence to the known analytic value, previously discussed 
m Ref. [3. We will return to a discussion of the analytic 
result later. Correspondingly, the real part gives shot 
noise, which is plotted in Fig. [2j and can be compared to 
the results in Ref. 0. Despite the rather larger transient 
behavior, the asymptotic behavior is quite clear, and is 
in perfect agreement with the previously discussed ana- 
lytic result for noise (which should be corrected for finite 
size effects [3j). We note that within the present method, 
the extraction of the (zero-frequency) noise is a lot easier 
than obtaining the correlation function directly, as noise 
is obtained directly from the large time limit of F(x) 
without having to perform a Fourier transformation in- 
cluding the associated finite size extrapolation [15j . 

Having examined a few time scans, we now turn to the 
FCS - but before examining the numerical plots let us 
give the analytical result for the CGF of this model. It 
has been long known [lj| that the IRLM is integrable - 
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that is all eigenvalues and eigenstates of the interacting 
many-body problem may be constructed via the Bethe 
ansatz. However it was only much later realized 17] 
that at the special interaction value given by the self- 
dual point, the voltage bias operator is also diagonal in 
this Bethe basis, allowing the exact solution to be found 
out of equilibrium via the methods of Ref . [IH by mapping 
to a boundary sine-Gordon (BSG) model. We note that 
this is not a generic feature of integrable models, a more 
general extension of the Bethe Ansatz to non-equilibrium 
systems is currently being discussed in the literature [l9| . 

The current and shot noise in the self-dual IRLM 
have previously been found and compared to numerical 
work in Refs. 0, S One can extend the previous analy- 
sis 2(J 21 1 to obtain all cumulants of charge transfer in 



the zero temperature limit. We find that the results are 
best expressed as a recurrence relation [22| 



r l (r v dGn 



(6) 



Taking the known expression for the current C\ [2| and 
resumming the cumulants, we find for small Vsd < V c : 
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and in the opposite circumstance Vsd > V c 
2oi(m) /Vsd 



F(x) =VsuJ2 



m>0 
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Here V c = -^zT' B is the convergence radius of each of 
the series, and T' B ~ (£') 4 / 3 is the natural energy scale 
of the problem. The non-universal parameter of propor- 
tionality relates the regularization of the field theory (on 
which the Bethe ansatz solution is based) to that of the 
original lattice model, and can be taken from previous 
works ■ This leaves zero free fitting parameters for all 
comparisons between analytic and numerical results in 
the current study. The series coefficients are given by 



a K (m) 



(-l) m+1 r (1 + Km) 
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We note here that the above results could be conjectured 
from only the self-duality of the model 10] . We also point 
out that the formal mapping between the IRLM and the 
BSG model is only proven at the self-dual point of the 
former model, which forces K = 4, 1/4 above. However 
a number of other models map onto the BSG model over 
a wider parameter range, to which the present results 
may be easily adapted. These include tunneling between 
quantum-Hall edges (20j . impurities in Luttinger liquids 
22] and I-V characteristics of Josephson junctions [23j]. 
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FIG. 3: [Color online] Cumulant generating function F(x) as 
determined by the DMRG, for various bias voltages Vsd < Vc. 
The data above is for t' = 0.3 for which V c = 0.64. The 
corresponding analytic results for the imaginary part is also 
shown in the plot. The analytic real part is too small to see 
on the current scale - the signal seen by the DMRG is almost 
entirely due to finite size effects. 



We now analyze the results, first considering the low 
bias regime. We plot the numerical data in Fig. [3] along 
with the analytic result fl7a|) . While the numerical data 
and analytic result show reasonable agreement in the 
imaginary parts, the real part suffers much more seriously 
from finite size effects, as was already noticed in Ref. y. 
In fact, the true values of the real part of F(x) hi the 
thermodynamic long time limit are very small, so almost 
all of the numerical signal is due to the finite size and 
finite measuring time. We note that while such large ef- 
fects might be a nuisance for numerical simulations, they 
are also extremely relevant for experiments where trans- 
port measurements through nano-devices may exhibit a 
similar behavior. By building a theory of such finite size 
effects and subtracting this from the data, some prelimi- 
nary agreement with the analytic result may be seen, in 
particular the ir periodicity that will be discussed shortly. 
Details of this will be presented elsewhere. 

Looking at the analytic expansion, Eq. I7a| we see 
that the backscattering current may be thought of as 
a sum of Poissonian processes of effective quas i-particlcs 

[iJ 111. As a re- 



of charges 2me where m is an integer 
suit of quantum interference, the effective 'probabilities' 
in this equivalent Poissonian process are not all positive, 
however from the periodicity of F(x), it is clear that the 
fundamental quasi-particles being scattered carry charge 
2e. This behavior was already hinted at by the Fano 
factor [3j, knowledge of the full CGF confirms this. 

We now turn to the large voltage regime: in Fig. [4] we 
plot numerical results along with the appropriate analytic 
curves from Eq. I7bl Again, we find very good agreement 
between the two results. As in the low bias voltage case, 
we can interpret the expansion (I7b|) as effective Poisso- 
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FIG. 4: [Color online] Comparison between numerical and 
analytic results for the cumulant generating function at bias 
voltages Vsd > V c . Data for different values of t' such that 
Vsu/T' B remains constant collapse onto the same curve, and 
show excellent agreement with the analytic result. 



nian processes, but this time with fractionally charged 
quasi-particles me/2 \2£ 



22], |24J, again confirming be- 



havior already hinted at by the Fano factor in Ref. |3|. 

While the function (|7bl) is actually 4-7T periodic, it is not 
clear from the numerical technique of explicitly adding 
the counting field that we can obtain a result that is 
not 2ir periodic. Furthermore, the numerical data is very 
messy for \ > tt making this parameter regime at present 
unobtainable to us. As we discuss shortly, the 47r period- 
icity formally arises due to F(x) becoming double valued 
- the branch may be chosen such that the function is 2tt 
periodic with discontinuities, or smoothly continued to be 
47T periodic. In an experimental situation, which branch 
is measured would depend on details of the setup 25] . 

Having found effective quasi-particles with charge e/2 
in the high bias limit, and with charge 2e in the low bias 
limit, we now discuss how to get from one to the other. 
In fact, this is easiest to see by studying the counting 
current I(x) — idF/dx, where we can actually resum 
the series and obtain the answer in closed form 



I{X) = ^sF 2 



{--1} {- -}-(v^ 2 



6' 6 J 



, (9) 



where V — Vsx>/V c is the normalized voltage, and 3F2 
is a hypergeometric function. This function has branch 
cuts along the lines \ — ±7^/2 from V = 1 to V = 00. 
This shows us technically how the CGF transforms from 
being n to 47r periodic at a critical bias V — 1. While the 
CGF for |x| > tt/2 undergoes a sharp transition at this 
voltage, the curves at low x, including all the cumulants, 
are smooth and merely show a crossover in their behavior. 

In summary, we have presented a numerical technique 
based on time-dependent DMRG from which one can 
extract the finite-bias CGF of the FCS for a generic 
quantum impurity problem. We have demonstrated the 



method for the self-dual IRLM, for which we have also ob- 
tained results from Bethe ansatz calculations, with good 
agreement. From this, we provide a scenario of obtaining 
charge fractionalization in transport experiments, which 
we relate to the analytic properties of the CGF. 

The numerics was performed at the XC2 and HC3 of 
the Steinbuch Centre of Computing, Karlsruhe, within 
project RT-DMRG. 
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